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Abstract 

We present the results of a variational Monte Carlo calculation of the 
exchange-correlation energy for a spin-polarized two-dimensional electron gas 
in a perpendicular magnetic field. These energies are a necessary input to the 
recently developed current-density functional theory Landau-level mixing is 
included in a variational manner, which gives the energy at finite density at 
finite field, in contrast to previous approaches. Results are presented for the 
exchange-correlation energy and excited-state gap at v = 1/7, 1/5, 1/3, 1, and 
2. We parameterize the results as a function of r s and v in a form convenient 
for current-density functional calculations. 
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I. INTRODUCTION 



Two-dimensional electron systems have been a subject of great interest for some time 
now, as fascinating phenomena such as the fractional quantum Hall effect and the Wigner 
solid have been studied intensively. More recently, technological progress has allowed ex- 
perimenters to fabricate extremely small structures, quantum dots, wires, and so on, that 
function as external potentials imposed on the 2-D electron gas. An extension of density 
functional theory, current density functional theory (CDFT) P]J2], has been used to study 
these phenomena ^M. The use of CDFT in this context is attractive, but the calculations 
to date have been hampered by the lack of an accurate exchange-correlation energy. The 
ones in use have relied on an interpolation between results at zero magnetic field as a func- 
tion of density |{7|||, and at infinite magnetic field and infinite density ||; i.e. the lowest 
Landau level approximation. The CDFT uses derivatives of the exchange-correlation energy 
with respect to density and magnetic field, so these interpolations will reflect only zero- and 
infinite- magnetic field and density properties of the system. In addition, since the lowest 
Landau level approximation requires that the density n and magnetic field B go to infinity 
with the ratio y/n/B held constant, a degree of freedom in the theory is lost. 

In this paper, we report the results of variational quantum Monte Carlo calculations of 
ground-state energies at various integer and fractional magnetic filling factors v. Our results 
include an estimate of Landau-level mixing, which allows us to calculate these quantities at 
finite density and magnetic field. In addition, we find the quasielectron-quasihole gaps at 
fractional filling factors, including the effect of Landau level mixing here as well. Using both 
these quantities yields a picture of the true exchange-correlation potential near v at varying 
density, including the FQHE-induced cusp. 

The CDFT generalizes the usual density functional theory by including the coupling of 
orbital currents to a magnetic field. The functionals in the theory depend not only on the 
particle density n(f) but on the paramagnetic current density j p (r) as well. (The physical 
current density j(r) is obtained using the continuity equation.) Solving the CDFT Kohn- 
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Sham equations || for the single-particle orbitals ip(r) then lets us define n(f) and j p (r) 
self-consistently, 



nir 



EI^WI 2 (i) 



lif) = ^ E [V£ (r) V^(r) - ^(r-)V^(f)] . (2) 

i 

The Kohn-Sham equations involve two exchange-correlation potentials: a scalar poten- 
tial V xc (r) and a vector potential A xc (r). These potentials, in turn, are functionals of the 
exchange-correlation energy E xc as follows, 

on(r) 
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where we have defined the vorticity 



v(f) = V x ( . (5) 
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The fact that enters in only through the vorticity ■u(r) ensures that the potentials are 
gauge-invariant. In addition, a local current density approximation (LCDA) can be defined 
in terms of the density n{f) and the vorticity v(r): 

EJMf),m= f^n^M^m (6) 



As with the LDA, this local energy e xc is the exchange-correlation energy of a uniform 
state, in this case the FQHE liquid state. Because the magnetic field in this state is uniform, 
translational invariance requires that the physical current density j(r) must be uniform and 
therefore zero. Then 

jfa) = Jp(*0 - —Mrl = 0, (7) 
mc 

where A(r) is the physical vector potential and n is the density. Taking the curl of both 
sides yields 



V x ( 3 -^-) = —V x A(r) = —B. (8) 
\ n J mc mc 

The left-hand side of (§) is just the vorticity, so using the definition of the magnetic filling 
factor v = 27rnohc/eB, we have 

2nhn 

v(r) = . (9) 

mv 

The vorticity in the uniform liquid state is then a function of the density no and filling 
factor v. Then in order to use the LCDA, we need only find the exchange-correlation energy 
of the uniform system 

s xc (r s , v) = e xc (n(r),v(r)) (10) 
as a function of density (or equivalently ion-disk radius r s = (7m ) _1 ^ 2 ) and filling factor v. 



II. VARIATIONAL MONTE CARLO METHOD 

We calculate the exchange-correlation energy e xc using a variational quantum Monte 
Carlo method. The method uses the Metropolis algorithm to evaluate the energy of a 
given wavefunction at various values of a variational parameter a. Since we are interested in 
including Landau-level mixing in the wavefunction only at certain fixed values of z/, where the 
liquid state is incompressible, we would like to use wavefunctions with suitable correlations 
built in already. The parameter a will be used to vary the amount of Landau-level mixing 
in the wavefunction. We can do this if we use lowest Landau level (LLL) wavefunctions for 
the incompressible states which are well-known: i.e. the Laughlin states for fractional v or 
Slater determinants for integer v. We obtain our variational wavefunction %p a by multiplying 
the LLL wavefunction by a Jastrow factor J a which lifts the wavefunction partially out of 
the LLL: 

(11) 
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where Rij is the distance between the electrons i and j, and u{r) is the so-called pseudopo- 
tential. 

Now, since ip a = ipoJa, where ?/>o is the LLL wavefunction, varying a varies the amount of 
Landau-level mixing. (The Jastrow factor forces ip a out of the LLL since it is non- analytic.) 
When a = we recover the LLL wavefunction. 

The pseudopotential u(r) is chosen to reproduce the zero-point motion of the plasmon 
modes in two dimensions f|; then 

u(r) ~ -±= (12) 



at large r. This pseudopotential has a singularity at the origin. In order to minimize the 
kinetic energy of the wavefunction we should require u{r) — > 1/y/F and ^ — > — ~ (^) ^ 
as r — > 0, where F is another variational parameter. A pseudopotential satisfying these 
constraints is given by 



1 — cxp 



r r \ 
F _ 2FJ 



(13) 



In practice, we have found that choosing F ~ 0.5 such that ^ — > — | ^ = — 1 minimizes 
the energy at arbitrary a, and thus r s , quite well. In all the calculations that follow, we have 
used F = 0.5, and varied only the parameter a. A comparison of pseudopotentials (0) and 
(|l3l) is shown in Figure 0. 

Whether our variational wavefunction is the ground state or an excited state, we wish 
to find its total energy 

E a = (a\H\a) (14) 

1 2 

= — (ali^la) H (alVla) , 

rf r s 

where H is the Hamiltonian, K is the kinetic energy operator, V is the potential energy 
operator and \a) is the many-particle state. The ion-disk radius r s factors completely out of 
both kinetic and potential parts. These two energies, (a|i^|a) and (a|V|a), are computed 
for a range of a. The energy E(r s ) is found at any given r s by minimizing E a with respect 



to a. This allows us to find E(r s ) in a corresponding range of r s without recomputing E a 
for each new r s . 

The potential energy (a|V|a) is easy to calculate by Monte Carlo. We write 

(a\V\a) = Jdn\i; a \ 2 V(n) (15) 

where Q = {fii, • • • , Ojv} and Qi is the position of the ith electron on the sphere. The 
integral is taken over all flj. We sample from the probability density \ip a \ 2 and sum the 
potential energy 

nn) = Eir- ( 16 ) 

i<j n iJ 

The kinetic energy sum is also straightforward, but more computationally demanding. 
The integral becomes 

(a\K\a) = [dn\^ a \ 2 ^. (17) 
J ip a 

As ip a will sometimes involve a projection of an already complex wavefunction to the LLL, 
we will have to take pains to make this computation efficient. If the projection is necessary, 
the cost of computing (Kip a )/ip a is of order N 3 , the most computationally expensive part 
of the code. 

A. Single-particle states and operators 



Our calculations are done on the sphere [11] rather than in the plane, as there are no 



boundaries and the calculations are more efficient. Since the calculations involve Landau- 
level mixing, which makes the particle density finite, we choose the radius of the sphere 
R = \^N/2, where N is the number of particles, to keep the area per particle constant. 
The magnetic monopole at the center of the sphere has magnetic charge S. We use units of 
length a = (vmo) -1 / 2 and units of energy e 2 /eae, where a# = h 2 e/m*e 2 and r s = a/as- The 
position of an electron on the sphere is described in spinor coordinates x — ( u , v ), where 
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and (8, 0) is the electron position in spherical coordinates. 

The kinetic energy of a single particle on the sphere can be described in terms of its 



angular momentum [10,11 



K 



R 2 



If we choose the gauge 



A = — — cot i 
R 



we find that the single-particle states are the monopole harmonics |12 | 



V —At „.S-m.,S+mr>S+m,S—m, 



UU — VV) 
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(21) 



M: 



S,l,m 



21 + 1 {i-sy.(i + sy. 
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An (l-m)\(l + m)\ 

Here / > S is the Landau-level index, and P^ ,l3 (x) is a Jacobi polynomial. (A bar over a 
variable denotes conjugation.) If we set I = S we recover the LLL states given in |JTT[ . 



The angular momentum operators given in |L1] are valid only in the LLL. Operators valid 
for all Landau levels must include the conjugate variables u and v. The angular momentum 
operators then become 



d d 
L+ = u- v— 

ov ou 

du dv 
l(_ d_ d_ -d__-d_ 
2 \ du dv du dv 



(22) 



These operators acting on the single-particle states Ys,i, m behave in the expected fashion; 
that is, L + raises m by one, L_ lowers m by one, and L z has eigenvalue m: 



L>+Ys,i, m = \J(l - m)(l + m + 1) Y St i t m+i 
L-Y Sjhm = y/(l + m)(l-m + l) Vs,i,m-i 
L z Ys i m — mYs i m . 



(23) 



Another set of operators exists on the sphere when a magnetic monopole is introduced: 

d d 
ov ou 

_ 1 / d d _ d _ _ d \ 
2\du dv du dv J 

These operators obey the commutation relations for angular momentum (SU(2)) just as the 
L operators do; i.e. [Ti,Tj] = e^Tfe. They commute with the L operators: [L, T] = 0, where 
L = (L x , L y , L z ), T = (T x ,T y ,T z ). Rather than changing the angular momentum quantum 
number m, however, they change the monopole charge S, 



T + Y s ,i, m = J(l-S)(l + S+ 1) Y s+1 



Lm 



T-Y s ,i,m = y/(l + S)(l-S+ 1) Y s _ ltljm (25) 
T z Ysj >m = S Ygj 

Note that T z — Cl • L, where (l is a unit vector indicating a position on the sphere. 

Because the monopole harmonics Ys t i >m form a complete basis for a fixed S, any wave- 
function if) which describes a particle moving on the sphere in the presence of a magnetic 
monopole of charge S must satisfy 

T z if) = Sif). (26) 

Now we can calculate the kinetic energy K in ([19]) more efficiently using the T operators. 
We use the facts that T 2 = L 2 and [T+, T_] = 2T Z to get 

L 2 = T 2 = - (T_T+ + T + T_) + T 2 

= T_T + + T 2 + T z (27) 
= T_T+ + S 2 + S. 

Then the kinetic energy operator becomes 

T_T + + S 

K = R2 , (28) 
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and T_T + acting on a wavefunction if) is enough to give us the kinetic energy. A wavefunction 
entirely in the LLL has T_T + ip = 0. 

These two sets of operators are analogous to the inter-Landau level operators a, a* and 
intra-Landau level operators b, in the planar geometry. This is easily seen by comparing 
the effects of these operators on the single-particle states in the plane to the T T and L T 
operators acting on Ys^ m on the sphere. 

B. Integer quantum Hall states 

The integer quantum Hall effect (IQHE) wavefunctions are commonly written as a Slater 
determinant in which all single-particle states up to the u th Landau level are filled. Any 
excitation of this state must involve raising one or more particles up one Landau level, 
producing an energy gap in the spectrum. This energy gap causes the many-particle state 
to be incompressible and creates the quantum Hall effect. There is no Landau- level mixing 
in this wavefunction, since each Landau level is either completely filled or completely empty. 
This wavefunction is strictly valid only at infinite magnetic field and infinite density. We 
introduce Landau-level mixing, and thus find the dependence of the total energy on the 
density, by attaching a Jastrow factor 

i>« = DJ a , (29) 

where D is the Slater determinant. The Landau-level mixing costs some kinetic energy but 
the cost is more than offset by a gain in potential energy, as correlations in the Jastrow 
factor allow the electrons to stay farther apart. 

The ground state is a uniform liquid, but excited states will create regions of charge 
excess and charge deficit. In order to create a free electron and hole in the excited state, we 
move these regions of charge excess and deficit as far apart as possible. On the sphere, this 
is best done by keeping the hole near the north pole and the electron near the south pole. 
We do this by removing an electron from the Y s ^-i state and placing it in the Xs^+i^+i state 



9 



in the first empty Landau level. The remaining interaction between hole and electron is 
removed by subtracting from the potential energy the interaction between an electron fixed 
at the south pole and a hole fixed at the north pole. 

The Monte Carlo algorithm proceeds by moving a single particle, chosen at random, 
then calculating the ratio of the new probability density \ip' a \ 2 to the old \ip a \ 2 - (Primes on 
a symbol indicate that it has been computed after the move.) If the new configuration Q' is 
more likely than the old Q, then the move is always accepted; if it is less likely, the move is 
accepted with probability 

Because the Jastrow factor is a pair product, the computational cost of an update after a 
single move is of order N. At first glance, the cost of an update to the Slater determinant 
is also of order N, since we can write 

D' = Y / MK)D mn , (31) 

m 

where we have moved the n th particle, 4> m (Q' n ) is the m th single-particle state after the 
move, and D mn is the appropriate cofactor. If the move is rejected, we need do no more 
computation. However, if the move is accepted, we now have to update all the cofactors 
with index n. 



Fahy et al |fL3|j , and earlier Ceperley et al |14|] have given an algorithm which computes 
these cofactors at a cost of order iV 2 . They begin the calculation by computing the inverse 
of the transposed Slater matrix, which gives the cofactors divided by the determinant 

D,k = • (32) 
The ratio of the new determinant to the old is 

q = ^ = Y,MK)D m n. (33) 

m 

The cofactors Djk are never calculated explicitly - they are kept in memory as the elements 
of the matrix Djk- These numbers are updated according to 
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Dai. la k = n 

D' ]k =\ _ (34) 



Computing the potential energy is now straightforward using (|i6|). 

The kinetic energy is computed by applying Ya ^i^t t° an d using Q28|). Because the 
T^ 1 are first-order differential operators, both the product rule and the chain rule apply, and 
we find 

T-T~T + ib 1 
- ' " Y. T, T/iD.U 



i) a DJ a 



E 



T~J a T?D , TrTfD 



J a D D 



+ "H*— (35) 



Tj D Tj J ex T% Tj J a 
D T 7 

The action of on the determinant is easily found, 



D 

T - T+D 

n n 

D 



Y,T±<P m (ttn)D mn (36) 

TO 

Y,T-T+<j) m {n n )D mn (37) 



The second quantity (0), when summed over all n, is just the excess kinetic energy of the 
determinant D, a constant independent of the electron positions. The first quantity (|36[) 
can be also be computed quickly, since the single-particle states are the monopole harmonics 
Ys,i,m(Qn) ■> an d we have already computed the matrix elements D nm . 

The action of T^- on the Jastrow factor is found most conveniently by defining two 
rotationally-invariant quantities on the sphere, 

nj = UiVj - ViUj (38) 
Sij = UiUj + ViVj. (39) 

Here |r^| is proportional to the chord distance on the sphere between particles i and j, 
and is related to the particle's center of mass. Their property of rotational invariance is 
verified by applying an arbitrary rotation 
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u v 
-v u 



( \ 

UqU + VqV 



V) 



UqV — VqU 



(40) 



substituting for u and v in (^) and fl39|), and using the identity uu + vv = 1 
The operators have a particularly simple effect on and Sy, 





Tij 


Ti 


= 


T~s- 


Tij = 


J- 1 It] 


= Sij 


T + s-- 
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± l s l] 



T + s-- 



'j 



(41) 



Other identities are found using 
Jastrow factor is 

rp± J 



-Tij and Sji = s^. Then the action of 7^ on the 



J, 



T n exp 



-aJ2 c ( r ij r ij) 

i<j 



(42) 



where c'(x) = dc/dx, the prime on the sum indicates that the sum is over all i ^ n, and 
we have defined a new pseudopotential function c(r 2 ) = u(r) equivalent to the old with a 
different argument. Applying (f41~[ ) gives explicit expressions. The quantity (T~T+J a )/J a is 
evaluated in a similar fashion. 



C. Fractional quantum Hall states 

It is well-known that the FQHE occurs at certain fractional filling factors v. In contrast 
to the IQHE, which is a single-particle phenomenon, the FQHE is caused by interactions 
between the particles. At these fractional fillings the particles can avoid each other best by 
staying in high angular momentum states relative to each other. Excitations in this state 
must have a gap, since any excitation has to lower at least one particle's angular momentum 
with respect to the others, using a finite amount of energy. 

The excitations appear as quasihole-quasielectron pairs, each carrying a fractional charge. 
The energy required to create the pair, the gap energy E g , determines the magnitude of the 
cusp that appears in the ground-state energy with respect to u, since the states just off the 
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FQHE filling factors v must have a finite number of quasiparticles. The magnitude of the 
cusp is related to the gap |T5j by 



'dE 




dE 




dv 


1>=V+ 


dv 


V = V — - 



(43) 



In this paper we will consider only the principal states at v = l/(p+ 1), p even. The 
Laughlin wavefunctions ipL describing the states on the sphere JTOf are written ip L = L p+1 
where 

L = JKuiVj - ViUj). (44) 

i<j 

We need an explicit wavefunction for the quasihole-quasielectron state as well. In order to 
do so we turn to the composite fermion picture, where 

iIjcf = DL P . (45) 

Here the Slater determinant D is a full Landau level when describing the principal states. 
Since at this filling D is a Van der Monde determinant, D = L and the two forms are 
equivalent. 

In the composite fermion picture, however, we think of D as describing the particles, and 
the remaining factor L p as attaching p magnetic flux quanta to each particle. To complete 
the picture we must then project the new wavefunction onto the LLL 

tpCF = V [DL P ] , (46) 

where V is the projection operator. To excite a quasiparticle pair, we follow the IQHE 
example by moving an electron in D from the LLL to the second Landau level. We again 
want the quasiparticles as far apart as possible, so we put the quasihole on the north pole 
and the quasielectron on the south pole. As before, we do this by moving the electron 
in Y s ,s-s to ^5,5+1,5+1. We then subtract from the potential energy the interaction of 
fractionally-charged particles on opposite poles to model a free quasihole and quasielectron. 

Unfortunately, projecting the composite fermion wavefunction to the LLL is, to the 
authors' knowledge, an unsolved problem. It has been done for a small (up to 12) number 

13 



of particles only by explicitly decomposing the wavefunction in the Slater determinant basis 
and keeping only those terms in the LLL [[HJ. There is, however, a trick recently used by 
Bonesteel |T7|] which can project a wavefunction with a single electron, or perhaps two or 



three, in the second Landau level to the LLL. If the occupancy of the second Landau level 
is at most one, for example, ipCF can be decomposed as 

$CF = lp°CF + 4>CF, ( 47 ) 

where ip^p is entirely in the LLL, and tpQ F has exactly one electron in the second Landau 
level. Letting the excess kinetic energy operator X) n T~T+ ac ^ 011 i>CF> then 

n 

=► (E T-lt - A)^cf = -^° CF (48) 

n 

is all in the LLL. Here A = 2S + 2 is the separation between Landau levels, and we have 
shifted the zero of the kinetic energy to ignore zero-point motion. Remembering the defini- 
tion in (|5BD , and normalizing by the unprojected wavefunction, we find 

4°CF 1 



tpCF tPcF^ 



S"T~T + - A 

/ j n n 



[Dm 



P^T-L T+D , (A D 



where we have used the fact that Tj~L = and 



n 



A D D = J2T-T+D. (50) 

n 

We have already found (T+D)/D. Writing L = Yii<j r ij an d using ([41]), we find 

S = -E'^- (51) 

i 'ni 

In order to introduce Landau-level mixing into the wavefunction, we again attach a 
Jastrow factor 



ij a = V [DL P ] J a . 
14 



(52) 



The potential energy, as before, can be found by averaging ([16]) over all configurations 
generated in the Monte Carlo code by \ip a \ 2 - 

The kinetic energy for the ground state can be found in a straightforward manner, since 
the ground-state variational wavefunction is 



ip a — L p+1 J a . 



(53) 



To calculate the kinetic energy, we use 



T~ T T+ T T~T+ T 

(j) _|_ ] j " n °a _|_ ' „ ' „ 



(54) 



where we have used the fact that T+L = 0. We have previously shown how to compute each 
of these terms. 

Computing the kinetic energy for the quasihole-quasielectron state is more involved, due 
to the presence of the projection operator V. The complexity arises from the fact that our 
implementation of V involves a product of two T operators, and the kinetic energy uses 
another product of two T operators. We need to calculate the quantity 



Em T m T+ (ipcF J° 



^CfJu m ^CF Ja 4>CF m Ja 

The only term we have not yet found is (T~i/;q F )/i/;cf- We have 

tMf t-v[dd>] 



(55) 



CF 



DLp 



T" 

DLp m 



p T-L T+D (A D 
A„ L D V A 



DL P 



pip - 1) T-L v T-L T+D 



P_ 

A 



E 



D 

T~T~L T+D T-L T'T+D 



(56) 



L 



D 



D 



+ 



A 



D 



A 



T+D T-L 

— h 

D L 



The two terms not yet computed are (T m T n L)/L and (T m T+D)/D. 
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We can find (T m T n L)/L in a straightforward manner, by applying (f4T|) and (|5TD. First 
we note that 

where we assume that n^i and we use the identity s^-s^ + r^f^ = 1. Then we have 



E 2- Er=l+ (58) 



1 — 5 / <5 2 



„■ „■ j.™ \ ' mi / 



The first two sums have been computed earlier, and the sum in the last term enters only 
when n = m. 

The second term (T~T+D)/D is computed easily when m = n, 

T~T + D 



D 

When m/nwe use the identity (derived in the appendix 



J2T-T+UVn)D in . (59) 



T~T + D 

m n 

D 



'£T+<f> i (n n )D in j hrTMttjD^j - (60) 

where the <fii{Qj) are again the single-particle wavef unctions. The N 3 scaling of the time 
required to compute the kinetic energy of ipQ F arises from (0), since ( pD| ) has to be computed 
for all N 2 pairs (m, n) and each computation scales as N. 



III. IMPLEMENTATION OF MONTE CARLO METHOD 

The implementation of the Monte Carlo procedure described above requires some adjust- 
ment in practice. The pseudopotential u(r), ([13]), if put directly on the sphere, will cause 
wild fluctuations in the kinetic energy calculation. The cause is found in the cusp at r = 2R, 
introduced in the pseudopotential by the mapping to the sphere. In effect, a ^-function in 



16 



the second derivative of the trial wavefunction appears if du/dr ^ at r = 2R. We remove 
the cusp by making a linear combination of w(r)'s with different arguments, 



u{r) = u(r) + u(4R — r). 



(61) 



Then 



du 




du 




du 


dr 


r=2R 


dr 


r=2R 


dr 



(62) 



4R-r=2R 



du/dr = at r = 2R, and the cusp is removed. The Jastrow factor is now 

J n 



exp 



-a^uflryl) 

i<j 



(63) 



where we have used the definition (0). With the variational parameter a appearing only 
outside the pseudopotential, we can compute potential and kinetic energies for a range of a, 
then do the minimization with respect to a for any number of values for r s without repeating 
the Monte Carlo calculation. 

Calculating the IQHE ground-state potential and kinetic energies by the Monte Carlo 
procedure in section |11 ij| is straightforward, since the trial wavefunction is simple. The 
variance in the result is small and scales as a/TV. The energy gap to the first excited state 
is more problematic, however. If we move an electron up to the next Landau level, creating 
an electron-hole pair, we must evaluate the pair energy by subtracting the ground state 
total energy from the excited state total energy. The variance in the Monte Carlo results for 
either state is then of the same order as the energy gap, and the results of the calculation are 
obscured by the statistical noise. Increasing N only worsens the problem, since the energy 
gap remains approximately constant with N, but the variance in the energies gets larger. 

We circumvent this difficulty in the manner recommended by Ceperley in PD|| , by sam- 
pling from the ground-state wavefunction and computing the excited state energy and ground 
state energy simultaneously. The two energies are then correlated: an upward fluctuation 
in the ground state energy always correlates with a nearly equivalent upward fluctuation in 
the excited state energy, for example, and the variance of the difference is small. 
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The computation of the excited state energy is implemented in the following way: the 
ground state energy is calculated in the Monte Carlo code as 

2H1PI 



1>. 



^nd 
gnd 



E 



vi 1 



nd ' 



(64) 



since the sample configurations are drawn from ip^ 1 (Q). (As before, we define Q = 
{Qi, ■ • ■ where fli is the position of the i th electron.) We write the excited state 

energy as 

l2 Hr* 



(a\H\a) ex = JdQ\r 



= k J dVt 

= «£ 

n 

where the normalization k is no longer 1 



^gnd 



gnd 



^r d (O r 



dfl 



E 



c x (^ 



4> c o 



(65) 



^I 1 
2 



gnd 



(66) 



Calculating the ground state energy for the FQHE states is comparatively straightforward 
as well, using (|6~4|). The FQHE ground state energies in this paper were all done in this 
manner. Because of the complexity of the excited state wavefunction, however, it is difficult 
to follow the recipe in (|65l). Instead, we sample from the unprojected wavefunction 



DL p J n . 



We compute both excited state and ground state wavefunctions similarly to 
^,gnd = <p [n gnd L p ] J a = L p+1 J a and ^ x = V [D CX L P ] J a . For the ground state, 

2 H^ nd 



(67) 
Define 



(a\H\a] 



gnd 



dtt 



1% 



gnd 
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un 1 2 
a I 



gnd 
a 



E 



2 



nd 



(68) 



and for the excited state 

(a\H\a) cx = Jd£l\P e 



E 



2 



Hr 



(69) 



The normalizations /t gn d and K ex are calculated in the same way as (|6q). 

We checked the results of this calculation by comparing the ground state energies found 
using (|68|) to the ground state energies found using (0). The energies were identical to 
within the statistical variations. 

The initial electron configuration for each Monte Carlo run must be chosen carefully to 
avoid underflow of the computer's arithmetic. The probability of generating a configuration 
from the trial wavefunction identical to that made by placing electrons at random on the 
sphere is so low that the initial calculation of the Slater determinant D will invariably 
underflow. To work around this problem we use an initial configuration which keeps the 
particles as far apart as possible. The probability of drawing this configuration from the trial 
wavefunction is much higher than that of a random configuration, so the initial computation 
of D will not underflow. 



This configuration is generated using some code borrowed from Jon Leech [J8]. It begins 
with six electrons; two on the north and south poles and four evenly distributed around 
the equator, forming the vertices of an octahedron. Each face is a triangle, which can be 
subdivided into four smaller triangles, as in Figure |2|. Placing electrons at the vertices of 
these new triangles generates a configuration of 18 electrons. Repeating the process generates 
another configuration of 66 electrons, and so on. The electrons are, of course, not directly 
on the vertices of the new triangles, but lie on the surface of the sphere, on a line joining 
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each triangle vertex to the origin. A few electrons can be removed from this configuration, 
as needed, without underflowing the computer arithmetic. 



IV. RESULTS 

We have obtained exchange-correlation energies e xc and gap energies A c h at v = 1/7, 
1/5, 1/3, 1, and 2 as a function of r s , using a model |20| of the electron-electron interaction 



^) - (70) 

Setting A = recovers the Coulomb interaction. The softer interactions modeled by ( |70|) 
more closely approximate the interaction experienced by electrons in a GaAs-AlGaAs het- 
erojunction. We have calculated e xc and A eh for A = 0, 0.2, 0.5, and 1.0. 

Numerical fits for our computed exchange-correlation energies are given in Tables |-|V. 
The Monte Carlo results are modeled using an order (2,2) rational function fit, 

_ _ a + air s + a 2 r 2 s 

e «" i + birs + b2r 2 ■ V 1 ) 

This form ensures that the energy e xc remains finite as r s — > oo. The Monte Carlo code 
produced good results from r s = to r s = 50 in most cases, although when e xc changed 
rapidly, the code gave results only up to r s = 35. This was true only for the Coulomb 
interaction at v — 1 and 2, with fitting parameters given in Table |. A typical set of results 
is shown in Figure [| 

The uncertainty in the numbers from the Monte Carlo code are given in the tables 
by showing the standard deviation in the values of the last digit in the parameter a in 
parentheses. The errors in e xc for higher A are much smaller than for the Coulomb case, 
since for a given change in inter-electron spacing, the softer interactions produce a smaller 
(sometimes much smaller) change in the total energy of the system. The fluctuations in 
the instantaneous energy in the Monte Carlo code are then smaller, and the corresponding 
results more accurate. 



20 



Numerical fits for the gap energies A e h are given in Tables |V|- |VII]| . Here we also use a 



rational function fit, but since the uncertainty in the gap energies are relatively large we use 
only an order (1,1) fit, 

A eh = 5±f!!-. (72) 
1 + hr s 

Gap energies for v = 1/3 and 1/7 are shown in Figures [| and ||, respectively. 

These energies were calculated with an 18-particle system, rather than the 66-particle 
system used to calculate the exchange-correlation energies. Fluctuations in the larger system 
tend to overwhelm the gap energy, since only one particle (or composite particle) is involved 
in creating the excitation whose energy we are measuring. 

The quasihole-quasielectron (or in the case of integer z/, hole-electron) interaction between 
particles of opposite charge on the north and south poles is removed by subtracting the 
interaction of point charges of appropriate magnitude at the same places. The uncertainty 
in the results is given, as before, in the numbers in parentheses beside the values for ao in 
Tables [V]- |V We keep three digits in the fitting parameters other than ao to remove errors 
in the fit caused by truncation. 



V. SUMMARY 

We have calculated exchange- correlation energies and energy gaps between ground state 
and first excited state for quantum Hall states at various filling factors v using a model 
electron-electron interaction which is realistic for electrons in GaAs-AlGaAs heterostruc- 
tures. Including Landau-level mixing in the calculation allows us to present results at finite 
magnetic field and electron density. Our results are presented in the form of numerical fits 
to simple rational functions, convenient for use in current-density functional theory calcu- 
lations. Calculating energy gaps as well as ground state energies allows us to model the 
energy near the quantum Hall state realistically. This work has been supported by the ONR 
and the NSF. 
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APPENDIX A: TWO SINGLE-PARTICLE OPERATORS ACTING ON A 

SLATER DETERMINANT 



Let the two single-particle operators be Oj and Oj, where i and j are particle indices 
and i ^ j. We wish to calculate (OiOjD)/ D where D is a Slater determinant composed 
of single-particle states <f>i(Clj). In the following will be the i, j th cofactor of the Slater 
matrix. Transposing and taking the inverse of the Slater matrix, we find Dij = Dij/D. The 
first operator acting on D gives 

d j D = y £O j M^j)Dkj- (Al) 

k 

The new object OjD can be thought of as a new Slater determinant D' . The new Slater 
matrix is identical to the old Slater matrix with the exception of the j th column, which 
changes from <f> m (£lj) to Oj4> m (Qj), l<m<N. Given the old inverse matrix elements D^, 
we would like to find the new elements D[-. The Sherman- Morrison formula (equations 
(2.7.4) and (2.7.5) of 0) gives 

d'h = D*-rrr ( A2 ) 

1 -h Zj 

where we define (following the notation in [21] ) 

z k = Yl DrnkUm (A3) 
m 

w k = ^D km v m = D kj . 

m 

Substituting, we find that 



z. 



k k 

& ' D - 1. (A4) 



D 

Now let the second operator Oj act on the new Slater determinant 
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OiD' 1 



D D 

jy_ 
~ 15 



D 



(A5) 



where we have used (|A2[) and (|A3|). Recognizing that D'/D = (OjD)/D and using ([AT 
we rewrite (|A5|) as 



OiOjD _ OjD OiD 



D 



D D 



Here, using (|A3|) and the fact that D ki is an inverse matrix element 



= Oj(j) k {Vlj)D k i - 5, 

k 

but i 7^ j by assumption, so we obtain the result 
OiOjD 



(A6) 



(A7) 



(A8) 



APPENDIX B: SPLINE INTERPOLATION BETWEEN QHE STATES 

Interpolation between energies Ej of adjacent quantum Hall states in filling factor v can 
be conveniently accomplished using cubic splines. The energy cusp at each state is included 
in the interpolation by requiring that the discontinuity in derivative of energy E be given by 
(fl3|), rather than requiring that it be zero, as normally done. If the filling factor (energy) of 
the j th state is given as Vj {Ej), and the discontinuities in derivative are Aj, the equations 
defining the spline parameters E'L following (3.3.7) of [21], are 



6 



E U 



u j+i ~ v ±± e". + - Uj E" 



6 



E j+1 ~ E 3 E 3 ~ E 3-l 



A, 



(Bl) 



Solving this system of equations for the E'J parameters, and using them following [21| section 
3.3, the interpolation in energy between states at a given r s is easily found. Figure || is an 
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example, showing the cusps at r s = 0. Of course, this interpolation ignores the heirarchy 
states completely. 



24 



REFERENCES 

[1] G. Vignale and Mark Rasolt, Phys. Rev. Lett. 59, 2360 (1987). 
[2] G. Vignale and Mark Rasolt, Phys. Rev. B 37, 10685 (1988). 
[3] G. Vignale, Phys. Rev. B 47, 10105 (1993). 
[4] M. Ferconi and G. Vignale, Phys. Rev. B 50, 14772 (1994). 
[5] P. Skudlarski and G. Vignale, Phys. Rev. B 48, 8547 (1993). 
[6] Mark Rasolt and F. Perrot, Phys. Rev. Lett. 69, 2563 (1992). 
[7] B. Tanatar and D. M. Ceperley, Phys. Rev. B 39, 5005 (1989). 
[8] D. Ceperley, Phys. Rev. B 18 3126 (1978). 

[9] D. Levesque, J. J. Weis, and A. H. MacDonald, Phys. Rev. B 30, 1056 (1984). 

[10] F. D. M. Haldane, Phys. Rev. Lett. 51, 605 (1983). 

[11] G. Fano, F. Ortolani, and E. Colombo, Phys. Rev. B 34, 2670 (1986). 

[12] T. T. Wu and Chen Ning Yang, Nucl. Phys. B 107, 365 (1976). 

[13] S. Fahy, X. W. Wang, and Steven G. Louie, Phys. Rev. B 42, 3503 (1990). 

[14] D. Ceperley, G. V. Chester, and M. H. Kalos, Phys. Rev. B 16, 3081 (1977). 

[15] A. H. Macdonald and S. M. Girvin, Phys. Rev. B 33, 4414 (1986). 

[16] J.K. Jain, Phys. Rev. Lett. 63 (1989) 199; Phys. Rev. B 40, 8079 (1989); ibid. 41, 7653 
(1990). 

[17] N.E. Bonesteel, Phys. Rev. B 51, 9917 (1995). 

[18] personal communication, Jon Leech (leech@cs.unc.edu) (1989) . 

[19] K. Binder, 

25 



Ed. , Monte Carlo methods in statistical physics, 2nd Ed. (Springer-Verlag, New 
York, 1986). 

[20] F. C. Zhang and S. Das Sarma, Phys . Rev. B 33, 2903 (1986). 

[21] W. H. Press et al, Numerical Recipes in C: The Art of Scientific Computing, 2nd 
Ed., (Cambridge University Press, New York, 1992). 



26 



FIGURES 

FIG. 1. Pseudopotential u(r) for various a with cusp condition satisfied and not satisfied 
(|i~2|). The solid lines correspond to (|i~3| ) and the dashed lines to (12). 
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FIG. 2. Subdivision of a triangular octahedron face into smaller triangles as part of the initial- 
ization procedure for the Monte Carlo code. 
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FIG. 4. Quasiparticle gap A e h at v = 1/3 as a function of 
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FIG. 5. Quasiparticle gap A e h at v = 1/7 as a function of r s . Landau-level mixing continues 
to have an effect here, although the effect on total energy is miniscule. 
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FIG. 6. Spline interpolation between energies at v = 1/3, 1/5, and 1/7, showing cusps in energy 
at principal FQHE states. Here r s = 0. The dotted line is a smooth cubic spline interpolation 
between FQHE states. 
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TABLES 

TABLE I. Parameters for a (2,2) rational function fit to the correlation energy e xc for the 
Coulomb interaction (A = 0). 



V 


a 


ai 


h 


a 2 x 10 3 


b 2 x 10 3 


2 


-0.864(2) 


-0.3883 


0.4221 


-19.85 


18.97 


1 


-0.888(2) 


-0.7518 


0.8281 


-41.28 


39.52 


1/3 


-1.004(1) 


-0.01561 


0.01466 






1/5 


-1.0355(8) 


-0.007768 


0.007241 






1/7 


-1.0507(5) 


-0.004958 


0.004602 







TABLE II. Parameters for a (2,2) rational function fit to the correlation energy £ xc with A — 0.2. 



V 


a 


ai 


h 


a 2 x 10 3 


b 2 x 10 3 


2 


-1.7220(1) 


-0.2843 


0.4221 


-4.081 


2.327 


1 


-1.7322(1) 


-0.093231 


0.053603 


-0.21638 


0.12166 


1/3 


-1.75085(5) 


-0.00677 








1/5 


-1.75623(5) 


-8.2 x 10~ 6 








1/7 


-1.75891(5) 


-3.9 x 10~ 6 
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TABLE III. Parameters for a (1,1) rational function fit to the correlation energy e xc with 
A = 0.5. 



V 


a 


ai bi 


2 


-3.22068(1) 


-0.0460183 0.0142792 


1 


-3.22281(1) 


-4.458 x 10~ 6 


1/3 


-3.22492(1) 


-2.133 x 10~ 7 


1/5 


-3.22540(1) 


-5.4 x 10~ 8 


1/7 


-3.22562(1) 


-2.0 x 10~ 8 



TABLE IV. Parameters for a linear fit to the correlation energy e xc with A = 1.0. 



V 


a 


ai 


2 


-4.81968(1) 


-6.26 x 10~ 7 


1 


-4.82005(1) 


-8.45 x lO" 8 


1/3 


-4.82034(1) 


-2.73 x 10" 9 


1/5 


-4.82040(1) 


-6.5 x 10~ 10 


1/7 


-4.82042(1) 


-2.3 x lO" 10 



34 



TABLE V. Parameters for a (1,1) rational function fit to the gap energy A e h for the Coulomb 
interaction (A = 0.0). 



V 


a 


ai 


6i 


2 


0.43(4) 


0.047 


0.321 


1 


0.80(8) 


0.0239 


0.176 


1/3 


0.27(3) 


0.00317 


0.0674 


1/5 


0.076(8) 


-6.15 x 10~ 6 


0.0107 


1/7 


0.040(5) 


3.62 x 10~ 4 


0.0187 



TABLE VI. Parameters for a (1,1) rational function fit to the gap energy A e h with A = 0.2. 



V 


a 


ai 




2 


0.32(5) 


0.0213 


0.164 


1 


0.56(8) 


0.00735 


0.0711 


1/3 


0.13(2) 


0.000847 


0.0239 


1/5 


0.046(7) 


-0.000177 


0.0 


1/7 


0.025(4) 


0.000251 


0.0155 
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TABLE VII. Parameters for a (1,1) rational function fit to the gap energy A e h with A = 0.5. 



V 


a 


ai 


h 


2 


0.20(4) 


0.00568 


0.0472 


1 


0.30(6) 


0.00233 


0.0185 


1/3 


0.049(9) 


0.000545 


0.0152 


1/5 


0.018(4) 


-1.12 x lO" 5 


0.0 


1/7 


0.009(2) 


-6043 x 10~ 6 


0.0 



TABLE VIII. Parameters for a (1,1) rational function fit to the gap energy A c h with A = 1.0. 



V 


a 


ai 


bi 


2 


0.14(3) 


0.00248 


0.0196 


1 


0.17(4) 


0.00118 


0.00819 


1/3 


0.022(5) 


-3.7 x 10- 6 


0.0 


1/5 


0.008(2) 


-4.4 x 10~ 7 


0.0 


1/7 


0.004(1) 


-1.9 x 10~ 7 


0.0 
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